Why is the DNA Denaturation Transition First Order? 
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We study a model for the denaturation transition of DNA in which the molecules are considered 
as composed of a sequence of alternating bound segments and denaturated loops. We take into 
account the excluded-volume interactions between denaturated loops and the rest of the chain by 
exploiting recent results on scaling properties of polymer networks of arbitrary topology. The phase 
transition is found to be first order in d = 2 dimensions and above, in agreement with experiments 
and at variance with previous theoretical results, in which only excluded-volume interactions within 
denaturated loops were taken into account. Our results agree with recent numerical simulations. 

PACS numbers: 87.14.Gg, 05.70.Fh, 64.10+h,63.70.+h 



Thermal denaturation or melting of double-stranded 
DNA is the process by which the two strands unbind 
upon heating. The nature of this transition has been 
investigated for almost four decades Jl|-01. Experimen- 
tally, a sample containing molecules of a specific length 
and sequence is prepared. Then the fraction of bound 
base-pairs as a function of temperature, referred to as 
the melting curve, is measured through light absorption, 
typically at about 260 nm. For heterogeneous DNA, 
where the sequence contains both AT and GC pairs, the 
melting curve exhibits a multistep behavior consisting of 
plateaus with different sizes separated by sharp jumps. 
These jumps have been attributed to the unbinding of 
domains characterized by different frequencies of AT and 
GC pairs. The sharpness of the jumps suggests that the 
transition from bound to unbound is first order. 

The early theoretical models [^]]3), which we refer 
to as Poland-Scheraga (PS) type models, consider the 
DNA molecule as composed of an alternating sequence 
of bound and denaturated states (see, e.g., Fig. |l|). Typ- 
ically a bound state is energetically favored over an un- 
bound one, while a denaturated segment (loop) is en- 
tropically favored over a bound one. Within the PS 
type models the segments which compose the chain are 
assumed to be non-interacting with one another. This 
assumption considerably simplifies the theoretical treat- 
ment and enables one to calculate the resulting free en- 
ergy. In the past the entropy of the denaturated loops has 
been evaluated by modelling them either as ideal random 
walks H or as self-avoiding walks ||. It has been found 
that within this approach the denaturation transition of 
DNA is continuous both in two and three dimensions and 
it becomes first order only above four dimensions. It was 
suggested |J that taking into account the interaction be- 
tween loops would sharpen the transition. 

A version of the PS model which includes the excluded- 
volume interactions between the various segments of the 
chain has been recently introduced jHJ. This extension 
makes the model analytically intractable, but numerical 
studies in three dimensions indicate that these interac- 
tions indeed drive the transition first order. 

In this Letter we consider analytically the effects of 





FIG. 1. 
model. 



Schematic representation of the Poland-Scheraga 



excluded-volume interaction between the various seg- 
ments of the chain. Although we treat this interaction 
only in an approximate way, we are able to give some 
insight into the unbinding mechanism and on the nature 
of the transition. Our approach makes use of recent im- 
portant results on the entropy of self-avoiding polymer 
networks pj]Jl2[ | . We find that this interaction drives the 
transition first order in d = 2, 3 and 4 — e dimensions. 
To proceed we introduce and discuss the PS model. A 
scaling argument which takes into account the interaction 
between the segments of the DNA chain is then presented 
and used to study the transition. Finally, we comment 
on disorder and compare our results to those obtained 
from more recent models of DNA denaturation [^-[tJ . 

The PS model considers two strands, made of 
monomers, each representing one persistence length of 
a single strand (ordinar ily ~ 40 A Q). For simplic- 
ity we take boundary conditions where the monomers at 
the ends of the molecule are always bound. All other 
monomers on the chain can be either bound or unbound 
to a specific matching monomer on the second chain. 
The binding energy Eq < is taken to be the same 
for all matching monomers. A typical DNA configu- 
ration is shown in Fig. [l]. It is made of sequences of 
bound monomers separated by denaturated loops. The 
statistical weight of a bound sequence of length k is 
io k = exp(—kEo/T), where T is the temperature [[IJJ. 
On the other hand the statistical weight of a denaturated 
sequence of length k is given by the change in entropy 
due to the added configurations arising from a loop of 
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length 2k. For large k this has the general form As k /k c , 
where s is a non-universal constant and the exponent c is 
determined by the properties of the loop configurations. 
For simplicity, we set A = 1. The model is most easily 
studied within the grand canonical ensemble where the 
total chain length L is allowed to fluctuate. The grand 
canonical partition function, Z, is given by 



Z 
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Vo(z)U L (z) 
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(1) 



where G(M) is the canonical partition function of a chain 
of length M, z is the fugacity and 



Lu k z k , (2) 



fe=l 



fc=l 



Vo(z) = 1+^0) and £/ L (z) = l + C/(z). Equation (g) can 
be verified by expanding the partition function as a series 
in U(z)V(z). The factors Vq{z) and Ul(z) properly ac- 
count for the boundaries. The average chain length, (L), 
is set by choosing a fugacity such that (L) = d In Z/d In z. 
We wish to evaluate the order parameter 9, i.e., the frac- 
tion of bounded pairs, as a function of temperature, in 
the limit (L) — ► oo. The average number of bounded 
pairs in a chain is given by (m) = dhiZ/dhiLo, so that 



= lim 
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(L) 
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(3) 



Here z* is the value of the fugacity in the limit (L) — > oo. 
This is the lowest value of the fugacity for which the par- 
tition function (||) diverges, i.e., for which U(z*)V(z*) — 
1. Using V(z) = ljz/(1 — loz) one has 



U{z*) = l/(coz*) - 1. 



(4) 



It is clear that the nature of the denaturation transi- 
tion is determined by the dependence of z* on to. The 
function U(z) is independent of lo. It is finite for z < 1/s 
and diverges when z > 1/s. On the other hand the func- 
tion 1/V(z) increases continuously as lo decreases (corre- 
sponding to an increase in T). Thus, as the temperature 
increases from T = 0, the fraction of attached monomers 
decreases and z* increases. However if z* reaches the 
value z c = 1/s (so that l/V(l/s) > U{l/s)) any fur- 
ther increase of the temperature leaves z* unchanged 
so that 9 = 0. Therefore the transition takes place at 
z* = z c = 1/s. Its nature is determined by the behaviour 
of U(z) in the vicinity of z c . This is controlled in turn 
by the value of the exponent c. There are three regimes: 

1. For c < 1, U(z c ) diverges, so that z* is an analytic 
function of lo and no phase transition takes place. 

2. For 1 < c < 2, U(z c ) converges but its derivative 
U'(z) diverges at z* = z c . Thus 9 oc 8z*/8lo van- 
ishes continuously (0 ~ \T - T c |(2-c)/(c-i)) at thc 
transition and the transition is continuous. 



3. For c > 2, both U(z) and its derivative converge at 
z* = z c and the transition is first order. Here, in 
contrast to the continuous case, the average size of 
a denaturated loop is finite at the transition. 

The value of the exponent c can be obtained by enumer- 
ating random walks which return to the origin, so that 
c = dv. For ideal random walks this yields c = d/2: thus 
there is no transition at d < 2, a continuous transition 
for 2 < d < 4 and a first order transition only for d > 4 
||. On the other hand, for self-avoiding random walks 
the excluded volume interaction modifies the exponent 
to c = 3/2 for d = 2 and c « 1.764 for d = 3. The 
transition is thus sharper, but still continuous, in three 
dimensions 

In the above treatment the entropy of a loop of length 
21 was taken to be of the form s l /l c . In determining the 
exponent c the interaction of the loop with the rest of 
the chain has been ignored. In the following we argue 
that, when the interaction of the loop with the rest of 
the chain is taken into account, the above form for the 
loop entropy is still valid, but the exponent c is modified. 
The excluded volume interaction of the loop with the rest 
of the chain is expected to lower the entropy of the loop 
yielding a larger value of the exponent c. 

We evaluate the exponent c in this situation by consid- 
ering the number of configurations of a loop of length 21 
embedded in a chain of length 2L (see Fig. 2). In reality 
the chain itself is composed of an alternating sequence of 
loops and bound segments of various sizes. However, in 
the present analysis we ignore the detailed structure of 
the rest of the chain and consider it as being in a bound 
state. The fact that the rest of the chain is made of both 
bound and unbound segments increases the excluded vol- 
ume interaction and is thus expected to increase further 
the value of c. This point will be discussed later. 

In order to evaluate the exponent c for our configura- 
tion, we use the results obtained by Duplantier JTT 12 
for the entropy of general polymer networks: Consider 
a polymer network composed of segments and vertices 
with an arbitrary given topology. The results state that 
for a network composed of segments with an average 
length La the number of configurations is proportional 
to oc s La L a g ~ 1 where 7q is an exponent which depends 
only on the topology of the network. It is given by 



7g = 1 — dvL 



N>1 



n N cr N , 



(5) 



where C is the number of independent loops, is the 
number of vertices with N outgoing legs, and o~n is an 
exponent related to such a vertex. The exponents ctjv are 
exactly known in d = 2 and in d = 4 — e to 0(e 2 ). 

Consider now the relevant topology to the problem dis- 
cussed in this work (Fig. 2). We are interested in enu- 
merating the configurations of the network in the limit 
l/L -C 1, when the loop size is much smaller than the 
rest of the chain. To do this we note that when only two 
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VI. ^QV3 , V1 

FIG. 2. The topology of the loop embedded in a chain. 
The distance along the chain from a vertex of type VI to the 
nearest vertex of type V3 is given by L. The distance between 
the vertices of type V3 is equal to I. 



. V4 V4 

FIG. 3. An extreme topology where the loop of length 21 is 
between two denaturated loops of size 2L each. The vertices 
of order 4 are denoted by VA. 



length scales determine the geometry of the network, as 
in the case depicted in Fig . 2, the number of configura- 
tions can be written as fll2fl 



c = 2 



s L+ \L + iy^- l g{l/L), 



(6) 



for large L and I. Here g(x) is a scaling function and 
7ioo P can be evaluated using eq. (||). For the topology 
considered above of a loop embedded in two segments 
(Fig. 2) we have: one loop, C = 1; two vertices of order 
one, n\ = 2, corresponding to the two free ends of the 
chain (denoted by VI in the figure); two vertices of order 
three, n 3 = 2 (denoted in the figure by V3). Using eq. (||) 
we obtain 



7ioo P = 1 - dv + 2<7i + 2er 3 . 



(7) 



We are interested in a loop size much smaller than 
the length of the chain, l/L <C 1. Clearly, in the limit 
l/L — > 0, the number of configurations should reduce to 
that of a single self-avoiding open chain, which, to leading 
order in L, is given by s L L~ t ~ 1 , where 7=1 + 2a\. This 
implies that, in the limit x <C 1, 



9(x) 



(8) 



Thus the Z-dependence of T, which gives the change in 
the number of configurations available to the loop, is 



(9) 



It is therefore evident that for large I and L and in the 
limit 1/1 <C 1 the partition sum is decomposed into a 
product of the partition sums of the loop and that of the 
rest of the chain |l5). The excluded volume effect be- 
tween these two parts reflects itself in the value of the ef- 
fective exponent c associated with the loop entropy. This 
result is very helpful since it enables one to extend the 
PS approach to the case of interacting loops. From Eq. 
^ one sees that the appropriate effective exponent c is 

c = 7 _ 7loop = dv - 2cr 3 . (10a) 

In d = 2 cr 3 = -29/64 Q and v = 3/4 yielding 

c = 2 + 13/32. (10b) 

In d = 4 - e to 0(e 2 ) one has 03 = -3e/16 + 9e 2 /512 
and v = 1/2(1 + e/8 + 15/4(e/8) 2 ) yielding 



5e 2 /256 



(10c) 



In d = 3 one may use Pade and Pade-Borel approxi- 
mations to obtain (73 « 0.175 12 which with the value 
v w 0.588 @ yielding c w 2.115. 

Equations ( |Io|) are our main result. The fact that the 
exponent c is found to be larger than 2 in d = 2,3 and 
4 — e indicates that the transition is first order in d > 2. 

As noted above this approach does not take into ac- 
count all the excluded volume interactions in the system 
since the rest of the chain is assumed to be in a bound 
state. In practice however the rest of the chain is by itself 
composed of loops and bound segments of various sizes. 
The existence of loops is expected to enhance the ex- 
cluded volume interaction. Taking into account all these 
configurations in the partition sum is a formidable task. 
However, we estimate below the effect of these interac- 
tions on the partition sum and show that it does not 
modify the main conclusion of this analysis, namely that 
the transition is first order. 

For a given loop-bound configuration of the rest of the 
chain it is straightforward to check that the I dependence 
of the partition sum (|9|) in the scaling limit is valid as long 
as the two bound segments which are attached to the loop 
under consideration are long. On the other hand when 
at least one of these segments is short and the adjacent 
open loop is large the effective exponent c is modified. In 
order to estimate this exponent we consider the extreme 
case where the rest of the chain is composed of two loops 
each of size 2L (see Fig. [|). A similar analysis to that 
presented above yields for the value of c, 



c = dv — 04 , 
= 2 + 11/16 



= 2 + e/4- 15e 2 /128 



(11a) 

in d = 2, (lib) 
in d = 4 - e, (11c) 



where the values (74 = —19/16 in d = 2 and (74 = 
— e/2 + ll(e 2 /8) 2 in d = 4 — e dimensions are used 
along with those of v. Using (74 w 0.46 obtained by Pade 
and Pade-Borel approximations gives in d = 3 the value 
c w 2.22. It is easy to verify that this scaling behavior 
is unchanged as long as one of the bound states con- 
nected to the loop under consideration is small (so there 
is at least one vertex of order four). Evidently c for this 
configuration is larger than that of the configuration in 
Fig. 2. One can show that for arbitrary loop-segment 
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configuration of the chain the effective c is given by ei- 
ther Eqs. ( |l0| ) or Eqs. (|ll]) depending on the length of 
the bound segments attached to the loop. Since in both 
cases the exponent c is found to be larger than 2 in d = 2, 
d = 4 — e and in d — 3, this analysis strongly suggests 
that the transition is first order for d > 2. The analysis 
presented above is valid as long as the length I of an indi- 
vidual loop is much smaller than the total length of the 
chain. The fact that the transition is found to be first 
order, so that the average loop size remains finite at the 
transition, makes the analysis self consistent. 

The model treated above assumes, for simplicity, that 
the monomers at the ends of the two chains are bounded. 
In reality the ends need not be attached and one should 
consider open boundary conditions. Similar analysis for 
the case where one of the ends of the chain is open yields 
the same result. Note that one end of the chains must 
remain attached for a bound state to exist. In the above 
approach entangled configurations with non equivalent 
knot topologies, which are possible in d — 3, are counted 
as they should when one end of the chains is open. 

Finally, we compare the model to recently studied 
models for the denaturation transition. These have fo- 
cused on the stacking energy generating a varying stiff- 
ness in the bound and denaturated state J0,||,0 or on 
the helical structure of the DNA molecule M . All these 
studies are within a directed polymer approach and thus 
ignore the effects studied in this Letter. 

In summary, we have shown that the reduction in the 
number of available configurations of denaturated loops 
due to the excluded-volume interactions between a de- 
naturated loop and the rest of the chain is sufficient to 
drive the denaturation transition first order. This is done 
by extending the PS type models to include these in- 
teractions. Our results have a direct implication on re- 
cent numerical simulations of the model [ 10| . The model 
studied here treats the binding energy between matching 
monomers is taken to be equal, although each base can 
bind only to its corresponding base on the other strand. 
This last assumption is reasonable when the heterogene- 
ity of the chain is taken into account. For short segments 
to bind the detailed sequences of the two segments have 
to match and this is less likely to take place in a hetero- 
geneous system. It is therefore important to consider the 
effect of the heterogeneity of the binding energy between 
the matching base pairs on the order of the transition. 
Also, another effect which has recently been emphasized 
is the varying stiffness of the phosphate backbone. The 
model considered above does not treat these effects. 

It is instructive to note that although the average loop 
size is finite at the transition, the loop size distribution is 
rather broad and decreases as l/l c for large I. For c < 3 
this yields a diverging variance of the loop size. It would 
be of great interest to test this prediction experimentally. 

We are indebt to M. J. E. Richardson and H. Orland for 
many inspiring discussions and their involvement earlier 
stages of our study of the DNA denaturation transition. 



We also thank B. Duplantier for very helpful discussions 
and M. E. Fisher for useful comments on the manuscript. 
LP wishes to pay a tribute to Raffaele (Peo) Tecce, who 
introduced him to the beauties of biological systems. The 
support of the Israeli Science Foundation is gratefully 
acknowledged. 



[1] 
[2] 

[3] 

[4] 

[5] 
[6] 



Michael Visiting Professor. On leave of absence from Di- 
partimento di Scienze Fisiche and Unita INFM, Univer- 
sity "Federico II", Complesso Monte S. Angelo, 1-80126 
Napoli (Italy). 

For a review see R. M. Wartell and A. S. Benight, Phys. 
Rep. 126, 67 (1985); O. Gotoh, Adv. Biophys. 16, 1 
(1983). 

For a review see D. Poland and H. A. Scheraga (eds.), 
Theory of Helix- Coil Transitions in Biopolymers (Aca- 
demic, New York, 1970); F. W. Wiegel, in C. Domb 
and J. L. Lebowitz (eds.): Phase Transitions and Critical 
Phenomena Vol. 7, Pg. 101 (Academic, New York, 1983). 
T. L. Hill, J. Chem. Phys. 30, 383 (1959); B. H. Zimm, J. 
Chem. Phys. 33, 1349 (1960); S. Lifson, J. Chem. Phys. 
40, 3705 (1964); M. Ya. Azbel, Phys. Rev. A 20, 1671 
(1979). 

M. Peyrard and A. R. Bishop, Phys. Rev. Lett. 62, 2755 
(1989). 

D. Cule and T. Hwa, Phys. Rev. Lett. 79, 2375 (1997). 
S. Cocco and R. Monasson, Phys. Rev. Lett. 84, 5178 
(1999). 



[71 N. Theodorak opoulos, T. Dauxois ans M. Peyrard, cond 



mat/0004487 



D. Poland and H. A. Scheraga, J. Chem. Phys. 45, 1456 
(1966); J. Chem. Phys. 45, 1464 (1966). 
[9] M. E. Fisher, J. Chem. Phys. 45, 1469 (1966). 
[10] M. S. Causo, B. Coluzzi and P. Grassberger, Phys. Rev. 

E 62, 3958 (2000). 
[11] B. Duplantier, Phys. Rev. Lett. 57, 941 (1986); J. Stat. 

Phys. 54, 581 (1989). 
[12] L. Schafer, C. von Ferber, U. Lehr and B. Duplantier, 

Nucl. Phys. B 374, 473 (1992). 
[13] M. T. Record, S. J. Mazur, P. Melanon, J. H. Roe, S. 
L. Shaner and L. Unger, Ann. Rev. Biochem. 50, 997 
(1981). 

[14] For simplicity, the configurational entropy related to the 
embedding of the chain of bound monomers in ambient 
space is neglected. One may check that it will not modify 
our main results. 

[15] This scaling analysis is also valid for chain configura- 
tions considered by J. des Cloizeaux (J. Physique 41, 
223 (1980)), where the probability of contact between 
two points along a self-avoiding chain is calculated. These 
configurations were also considered in jnj. It should be 
noted that the simple approach described above holds as 
long as, in the limit I — > 0, no higher-order vertices (in 
the sense discussed above) are generated. 



4 



